Metadata¶
Preprocess cranberry bogs to remove unwanted artifacts and aggregate by proximity
Author: Adrian Wiegman
Date Modified: 10/31/2023
Creates a layer of cranberry bog polygons from raw input data
Main Program¶
# iphython options
# delete variables in workspace
%reset -f
#places plots inline
%matplotlib inline
#automatically reloads modules if they are changed
%load_ext autoreload
%autoreload 2
# this codeblock sets up the environment from jupyter notebooks
setup_notebook = "C:/Users/Adrian.Wiegman/Documents/GitHub/Wiegman_USDA_ARS/Cran_Q_C/2_gis/scripts/_Setup.ipynb"
%run $setup_notebook # magic command to run the notebook
The autoreload extension is already loaded. To reload it, use: %reload_ext autoreload *** loading python modules... `module_list` contains names of all loaded modules ...module loading complete *** loading user defined functions... type `fn_`+TAB to for autocomplete suggestions the object `def_list` contains user defined function names: fn_get_info fn_arcgis_table_to_df fn_arcgis_table_to_np_to_pd_df fn_run_script_w_propy_bat
Signature: fn_get_info(name='fn_get_info') Docstring: returns the source information about a given function name Source: def fn_classify_wetlands(x): """ Consolidates DEP wetlands into fewer categories """ if "MARSH" in x: return("MARSH") elif "WOOD" in x: return("FORESTED SWAMP") elif "SHRUB" in x: return("SHRUB SWAMP") elif len(x) == 0: return("NON WETLAND") else: return(x) File: Dynamically generated function. No source code available. Type: function
fn_try_mkdir
fn_hello
fn_recursive_glob_search
fn_regex_search_replace
fn_regex_search_0
fn_arcpy_table_to_excel
fn_agg_sum_df_on_group
fn_add_prefix_suffix_to_selected_cols
fn_calc_pct_cover_within_groups
fn_buildWhereClauseFromList
fn_FA_to_Q
fn_alter_field_double
fn_return_float
fn_classify_wetlands
use ??{insert fn name} to inspect
for example running `??fn_get_info` returns:
you can also use `fn_get_info(name={insert fn name})` to do the same thing as `??{insert fn name}`
...function loading complete
***
managing project directories:
working directory (`wdr`): C:\Workspace\Geodata\Cran_Q_C\
temp dir (`tdr.name`): C:\Users\ADRIAN~1.WIE\AppData\Local\Temp\ArcGISProTemp29724\tmpac5oidy7
input file directory (`idr`): C:\Workspace\Geodata\Cran_Q_C\inputs
output dir (`odr`): C:\Workspace\Geodata\Cran_Q_C\outputs
***
setting up arcpy environment...
default geodatabase path (`gdb_default_path`): C:\Workspace\Geodata\Cran_Q_C\Default.gdb
active map (`mx`) name: Map2_prep_gw
output coordinate system: NAD_1983_UTM_Zone_19N
... env setup complete
+++SETUP COMPLETE+++
#wdr = r"C:\Workspace\Geodata\Cran_Q_C"
print(wdr)
C:\Workspace\Geodata\Cran_Q_C\
# Select file geodatabase to store results\n",
#gdb = "Cran_Q_C.gdb"
#ap.env.workspace = os.path.join(wdr,gdb)
# copy raw bog data to the default geodatabase in working directory
in_feat = r"C:\Workspace\Geodata\Cran_Q_C\inputs\WMAbogsDRAFT2013\WMAbogsDRAFT2013.shp"
bogs_raw = "bogs_raw"
arcpy.management.CopyFeatures(
in_features=in_feat,
out_feature_class=bogs_raw)
Messages
Succeeded at Wednesday, December 4, 2024 10:04:58 AM (Elapsed Time: 1.14 seconds)
process cranberry bogs to split features.¶
Example of unwanted artifact of cranberry bog layer.
Polygons that have the same owner are joined together by thin corridors. This can be problematic if the bog overlaps with a stream reach (denoted by the green point). These corridors sections must be removed. The corridors can be removed by putting a negative buffer on the the polygons. Additionally it may be helpful to split polygons that are far apart. And to combine polygons that are close together.
# Negative buffer around perimeter to remove small corridors (this removes links bewteen bogs with the same permit that are not hydrologically connected)
buff_dist = 3
bogs_neg = "bogs_buff_neg{}m".format(buff_dist)
arcpy.analysis.PairwiseBuffer(
in_features=bogs_raw,
out_feature_class=bogs_neg,
buffer_distance_or_field="-{} Meters".format(buff_dist),
dissolve_option="NONE",
dissolve_field=None,
method="PLANAR",
max_deviation="0 Meters"
)
--------------------------------------------------------------------------- ExecuteError Traceback (most recent call last) Cell In[29], line 4 2 buff_dist = 3 3 bogs_neg = "bogs_buff_neg{}m".format(buff_dist) ----> 4 arcpy.analysis.PairwiseBuffer( 5 in_features=bogs_raw, 6 out_feature_class=bogs_neg, 7 buffer_distance_or_field="-{} Meters".format(buff_dist), 8 dissolve_option="NONE", 9 dissolve_field=None, 10 method="PLANAR", 11 max_deviation="0 Meters" 12 ) File C:\Program Files\ArcGIS\Pro\Resources\ArcPy\arcpy\analysis.py:1298, in PairwiseBuffer(in_features, out_feature_class, buffer_distance_or_field, dissolve_option, dissolve_field, method, max_deviation) 1296 return retval 1297 except Exception as e: -> 1298 raise e File C:\Program Files\ArcGIS\Pro\Resources\ArcPy\arcpy\analysis.py:1281, in PairwiseBuffer(in_features, out_feature_class, buffer_distance_or_field, dissolve_option, dissolve_field, method, max_deviation) 1277 from arcpy.arcobjects.arcobjectconversion import convertArcObjectToPythonObject 1279 try: 1280 retval = convertArcObjectToPythonObject( -> 1281 gp.PairwiseBuffer_analysis( 1282 *gp_fixargs( 1283 ( 1284 in_features, 1285 out_feature_class, 1286 buffer_distance_or_field, 1287 dissolve_option, 1288 dissolve_field, 1289 method, 1290 max_deviation, 1291 ), 1292 True, 1293 ) 1294 ) 1295 ) 1296 return retval 1297 except Exception as e: File C:\Program Files\ArcGIS\Pro\Resources\ArcPy\arcpy\geoprocessing\_base.py:532, in Geoprocessor.__getattr__.<locals>.<lambda>(*args) 530 val = getattr(self._gp, attr) 531 if callable(val): --> 532 return lambda *args: val(*gp_fixargs(args, True)) 533 else: 534 return convertArcObjectToPythonObject(val) ExecuteError: Failed to execute. Parameters are not valid. ERROR 000840: The value is not a Feature Layer. WARNING 000725: Output Feature Class: Dataset C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_buff_neg3m already exists. Failed to execute (PairwiseBuffer).
# positive buffer to return to original area
bogs_neg_pos = bogs_neg+"_pos{}m".format(buff_dist)
arcpy.analysis.PairwiseBuffer(
in_features=bogs_neg,
out_feature_class=bogs_neg_pos,
buffer_distance_or_field="{} Meters".format(buff_dist),
dissolve_option="NONE",
dissolve_field=None,
method="PLANAR",
max_deviation="0 Meters"
)
# split all features containing multipart polygons (unconnected) into single part features (one feature per polygon).
bogs_split = bogs_neg_pos+"_split"
arcpy.management.MultipartToSinglepart(bogs_neg_pos, bogs_split)
# aggregate polygons of bogs within 30 meters and sharing the same FID.
# note that this modifies the shape of polygons by creating new corridors that link between adjacent polygons
# this increases total area to the of all polygons in the layer (which can lead to overestimates of potentially restorable bog area).
agg_dist = 30 # meters
bogs_agg = "bogs_agg"
arcpy.cartography.AggregatePolygons(
in_features=bogs_split,
out_feature_class=bogs_agg,
aggregation_distance="{} Meters".format(agg_dist),
minimum_area=None,
minimum_hole_size="0 SquareMeters",
orthogonality_option="NON_ORTHOGONAL",
barrier_features=None,
out_table="{}_Tbl".format(bogs_agg),
aggregate_field="ORIG_FID"
)
# use identity to get FID of aggregated bogs
bogs_ident = "bogs_ident"
arcpy.analysis.Identity(
in_features=bogs_split,
identity_features=bogs_agg,
out_feature_class=bogs_ident,
join_attributes="ONLY_FID",
cluster_tolerance=None,
relationship="NO_RELATIONSHIPS"
)
# dissolve all features with that have the same value for FID of bog_agg
# this groups feature within 30m but without adding area
bogs_ident_diss = "bogs_ident_dissolved"
arcpy.management.Dissolve(
in_features=bogs_ident,
out_feature_class=bogs_ident_diss,
dissolve_field="FID_{}".format(bogs_agg),
statistics_fields=None,
multi_part="MULTI_PART",
unsplit_lines="DISSOLVE_LINES",
concatenation_separator=""
)
# use identity to get original raw bog FIDs
bogs_ident_diss = "bogs_ident_dissolved"
bogs_raw = "bogs_raw"
bogs_agg_diss_ident = "bogs_agg_diss_ident"
arcpy.analysis.Identity(
in_features=bogs_ident_diss,
identity_features=bogs_raw,
out_feature_class=bogs_agg_diss_ident,
#join_attributes="ONLY_FID",
cluster_tolerance=None,
relationship="NO_RELATIONSHIPS")
# Add geometry attributes to table
arcpy.management.CalculateGeometryAttributes(
in_features="bogs_agg_diss_ident",
geometry_property="bogsurf_area_m2 AREA;bogsurf_perim_m2 PERIMETER_LENGTH",
length_unit="METERS",
area_unit="SQUARE_METERS",
coordinate_system=None,
coordinate_format="SAME_AS_INPUT")
# Remove unwanted artifacts of preprocessing.
# select and delete any polygons that do not overlap with original bogs
# these areas are noted by FID_bogs_raw = -1
bogs_agg_diss_ident = "bogs_agg_diss_ident"
arcpy.management.SelectLayerByAttribute(
in_layer_or_view=bogs_agg_diss_ident,
selection_type="NEW_SELECTION",
where_clause="FID_bogs_raw = -1",
invert_where_clause=None)
arcpy.management.DeleteRows(
in_rows=bogs_agg_diss_ident)
bogs_agg_diss_ident = "bogs_agg_diss_ident"
arcpy.management.SelectLayerByAttribute(
in_layer_or_view=bogs_agg_diss_ident,
selection_type="CLEAR_SELECTION")
# remove rows where area is <100m
arcpy.management.SelectLayerByAttribute(
in_layer_or_view=bogs_agg_diss_ident,
selection_type="NEW_SELECTION",
where_clause="Shape_Area < 100")
arcpy.management.DeleteRows(in_rows=bogs_agg_diss_ident)
pink areas show farm areas after splitting but before aggregation red areas with black outline show area after aggregation green points are the points of maximum flow accumulation within each aggregated cluster of beds
copy = arcpy.management.CopyFeatures(
in_features=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_raw",
out_feature_class=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_raw_copy",
config_keyword="",
spatial_grid_1=None,
spatial_grid_2=None,
spatial_grid_3=None)
arcpy.management.CalculateField(
in_table=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_raw_copy",
field="FID_bogs_raw",
expression="!OBJECTID!",
expression_type="PYTHON3",
code_block="",
field_type="",
enforce_domains="NO_ENFORCE_DOMAINS")
copy = arcpy.management.CopyFeatures(
in_features=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_agg",
out_feature_class=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_agg_raw",
config_keyword="",
spatial_grid_1=None,
spatial_grid_2=None,
spatial_grid_3=None)
arcpy.management.JoinField(
in_data="bogs_agg_raw",
in_field="ORIG_FID",
join_table=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_raw",
join_field="OBJECTID",
fields=None,
fm_option="NOT_USE_FM",
field_mapping=None,
index_join_fields="NO_INDEXES")
arcpy.analysis.SummarizeWithin(
in_polygons=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_agg_raw",
in_sum_features="bogs_agg_diss_ident",
out_feature_class=r"C:\Workspace\Geodata\Cran_Q_C\Default.gdb\bogs_agg_raw_SummarizeWithin",
keep_all_polygons="KEEP_ALL",
sum_fields="bogsurf_area_m2 Sum;bogsurf_perim_m2 Sum",
sum_shape="NO_SHAPE_SUM",
shape_unit="SQUAREKILOMETERS",
group_field=None,
add_min_maj="NO_MIN_MAJ",
add_group_percent="NO_PERCENT",
out_group_table=None)